# -*- coding: utf-8 -*-
"""
Created on Sat May 15 13:09:56 2021

@author: Apple
"""

import math
import numpy as np
import sympy as sp
import pandas as pd
import os
import sys

class Logger(object):
    def __init__(self, filename="Default.log"):
        self.terminal = sys.stdout
        self.log = open(filename, "a")
 
    def write(self, message):
        self.terminal.write(message)
        self.log.write(message)
 
    def flush(self):
        pass

def calvalue(y,x,x0): # 函数值求解

    fx = float(y.evalf(subs={x:x0}))
    return fx

def derfunction(y,x): # 输出导函数表达式及在x0处的值
    dy = sp.diff(y,x)
    return dy

def newton(y,x,x0,eps,m): # y为函数表达式,x为自变量符号,x0为初始点,eps为精度,m为重根个数 

    fx = calvalue(y,x,x0) # 代入y，求解出x0对应的函数值
    
    y1 = derfunction(y,x) # 求解出导函数表达式
    f1x = calvalue(y1,x,x0) # 代入y1，求解出x0对应的函数值
    
    x1 = x0-m*fx/f1x # 计算xk+1
    
    k = 0
    K = []
    # t=1   
    # while abs(calvalue(y,x,x1))>=abs(calvalue(y,x,x0)):# 牛顿下山法
    #     t = t/2
    #     x1 = x0 - t*fx/f1x
    
    while abs(x1-x0)>=eps:# 迭代
        K.append([k,x0,fx])
        k+=1
        x0 = x1
        fx = calvalue(y,x,x0)
        f1x = calvalue(y1,x,x0)
        x1 = x0-m*fx/f1x
        
    K.append([k,x0,fx])
    
    return k,x0,fx,K

def calMR(y,x,x0,a): # 计算重根个数 
    b = 0
    fx = calvalue(y,x,x0)
    # print(fx)
    while(abs(fx)<0.001 and b<a):
        b+=1
        y = derfunction(y,x)
        fx = calvalue(y,x,x0)
    return b
        

def main():
    
    type = sys.getfilesystemencoding()
    sys.stdout = Logger('scheme3.txt')
    
    x = sp.symbols('x')
    # y = 27*x**3+54*x**2+36*x+8 # 定义函数表达式
    # y = 36*x**4-12*x**3+37*x**2-12*x+1
    y = 2*sp.exp(x-1)-x**2-1
    # y = sp.log(3-x)+x-2
    a = 4
    
    print('函数表达式为: f(x) =',y)
    
    x0 = 0 # 设置初始值
    eps = 0.0000001 # 设置精度
    m = 1
    print(f'计算精度为{eps}')
    
    k,x0,fx,K = newton(y,x,x0,eps,m)

        
    print(f'经过{k}次迭代后，xk的值为{x0}，此时f(xk)={fx}')
    
    b = calMR(y,x,x0,a)
    print('重根个数为:',b)
    
    x0 = 1 # 重置初始值
    print('改进牛顿方法求得结果如下：')
    k,x0,fx,K = newton(y,x,x0,eps,b)
    data = pd.DataFrame(K,columns=["k", "xk", "f(xk)"])
    print(data)
    print('前向误差为：',abs(fx),'后向误差为：',abs(2-x0))
    
    
if __name__ == "__main__":
    main()
